library(tidyverse)
library(rnaturalearth)
library(sp)
library(eurostat)
library(haven)
library(modelsummary)
library(vtable)
library(labelled)
library(kableExtra)
library(gtsummary)
library(viridis)
library(hrbrthemes)
library(naniar)

# set the working directory
#setwd("")

#this would have been impossible without https://toscano84.github.io/2018/07/using-r-to-analyse-the-german-federal-election/
#although the data are slightly different, this was extremely helpful

data_env <-readRDS("data_env_ess8.rds")

table(data_env$regunit)

#for the moment, I'm going by Country


guide_legend <- guide_colorbar(
  direction = "horizontal",
  barheight = unit(2, units = "mm"),
  barwidth = unit(50, units = "mm"),
  draw.ulim = F,
  title.position = 'top',
  title.hjust = 0.5,
  label.hjust = 0.5)

#base map

geodata <- get_eurostat_geospatial(
  output_class = "sf",
  resolution = "3",
  nuts_level = 1,
  year = 2016
)

geodata <- geodata %>% 
  replace_with_na(list(NUTS_ID=c("FRY", "ES7", "PT2", "PT3"))) %>% 
  dplyr::filter(!is.na(NUTS_ID)) %>%
  glimpse()

ggplot(geodata) +
  geom_sf(  col = "black", lwd = 0.3 ) +
  xlab(NULL) + ylab(NULL) +
  ggtitle("") +
  theme(axis.text.x = element_blank(), #removes all the ggplot line ticks
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        rect = element_blank())


env_pop <- data_env %>% 
  dplyr::select(envpop, cntry) %>%
  group_by(cntry) %>%
  summarize(mean_envpop=mean(envpop, na.rm=T), N=n()) %>%
  rename(CNTR_CODE=cntry) %>% 
  glimpse()

env_pop$CNTR_CODE[env_pop$CNTR_CODE == "GB"] <- "UK"

summary(env_pop$mean_envpop)

#join map to ESS dat
geodata <- geodata %>% #base map
  left_join(env_pop, by="CNTR_CODE")

#merge county avg into base map


summary(geodata$mean_envpop)


geodata <- geodata %>% 
  dplyr::filter(!is.na(mean_envpop)) %>%
  dplyr::filter(!is.na(NUTS_ID)) %>%
  group_by(CNTR_CODE) %>%
  summarize(mean_envpop=mean(mean_envpop, na.rm=T), N=n()) %>%
  glimpse()

geodata %>% 
  ggplot()+
  geom_sf( aes(fill=mean_envpop), col="white", lwd = 0.1 )+
  #geom_sf_text(aes(label=NAME_LATN), size=.5, col="black")+
  xlab(NULL) + ylab(NULL) +
  ggtitle("") +
  theme(axis.text.x = element_blank(), #removes all the ggplot line ticks
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        rect = element_blank())+
  expand_limits(x=25)+
  scale_fill_gradient2(
    low = "#F2F2F2",
    mid= "#A4A4A4",
    high = "#000000",
    midpoint = 1.55,
    name = "Environmental Populism",
    guide = guide_legend)
ggsave("figures/eumap.png", bg ="white", width = 6.5, height = 4.33)





#Europe is being a pain to do subnationally, trying Germany

geodata1 <- get_eurostat_geospatial(
  output_class = "sf",
  resolution = "3",
  nuts_level = 1,
  year = 2016
)

DE <- geodata1 %>% dplyr::filter(CNTR_CODE=="DE")


env_pop_de <- data_env %>% 
  dplyr::filter(cntry=="DE") %>% 
  dplyr::select(envpop, region) %>%
  group_by(region) %>%
  summarize(mean_envpop=mean(envpop, na.rm=T), N=n()) %>%
  rename(NUTS_ID=region) %>% 
  glimpse()


DE <- DE %>% 
  left_join(env_pop_de, by="NUTS_ID")


DE %>%
  ggplot()+
  geom_sf( aes(fill=mean_envpop), col="white", lwd = 0.1 )+
  #geom_sf_label(aes(label=NAME_LATN), size=.5, col="black")+
  xlab(NULL) + ylab(NULL) +
  ggtitle("") +
  theme(axis.text.x = element_blank(), #removes all the ggplot line ticks
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        rect = element_blank())+
  labs(title = "")+
  scale_fill_gradient2(
    low = "#F2F2F2",
    mid= "#A4A4A4",
    high = "#000000",
    midpoint = 1.55,
    name = "Environmental Populism",
    guide = guide_legend)
ggsave("figures/envpop_de.png", bg="white", width = 6.5, height = 4.33)


